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PROTECTED POLYMORPHISMS AND EVOLUTIONARY 
STABILITY OF PATCH-SELECTION STRATEGIES IN STOCHASTIC 

ENVIRONMENTS 

STEVEN N. EVANS, ALEXANDRU HENING, AND SEBASTIAN J. SCHREIBER 

Abstract. We consider a population living in a patchy environment that varies 
stochastically in space and time. The population is composed of two morphs (that 
is, individuals of the same species with different genotypes). In terms of survival 
and reproductive success, the associated phenotypes differ only in their habitat selec- 
tion strategies. We compute invasion rates corresponding to the rates at which the 
abundance of an initially rare morph increases in the presence of the other morph es- 
tablished at equilibrium. If both morphs have positive invasion rates when rare, then 
there is an equilibrium distribution such that the two morphs coexist; that is, there is 
a protected polymorphism for habitat selection. Alternatively, if one morph has a neg- 
ative invasion rate when rare, then it is asymptotically displaced by the other morph 
under all initial conditions where both morphs are present. We refine the character- 
ization of an evolutionary stable strategy for habitat selection from [Schreiber, 2012] 
in a mathematically rigorous manner. We provide a necessary and sufficient condition 
for the existence of an ESS that uses all patches and determine when using a single 
patch is an ESS. We also provide an explicit formula for the ESS when there are two 
habitat types. We show that adding environmental stochasticity results in an ESS 
that, when compared to the ESS for the corresponding model without stochasticity, 
spends less time in patches with larger carrying capacities and possibly makes use of 
sink patches, thereby practicing a spatial form of bet hedging. 



1. Introduction 



Habitat selection by individuals impacts key attributes of a population including its 
spatial distribution, temporal fluctuations in its abundance, and its genetic composition. 
In environmentally heterogeneous landscapes, individuals selecting more favorable habi- 
tats are more likely to survive or reproduce. As population densities increase in these 
habitats, individuals may benefit by selecting previously unused habitats. Thus, both 
environmental conditions and density-dependent feedbacks generate selective pressures 
on habitat selection. Under equilibrium conditions, spatial heterogeneity can select for 
populations exhibiting an ideal-free distribution-equal per-capita growth rates in all oc- 
cupied patches and lower per-capita growth rates if individuals moved into unoccupied 
patches Fretwell and Lucas , 1969 . Under non-equilibrium conditions, spatial-temporal 



heterogeneity can select for individuals occupying sink habitats in which the per-capita 
growth rate is always negative Holt, 1997, Jansen and Yoshimura, 1998 . Environ- 



mental heterogeneity can also promote coexistence of genotypes only differing in their 
habitat choices [Jaenike and Holt 1991] . Despite significant advances in the mathemati- 
cal theory for habitat selection under equilibrium conditions, a mathematical theory for 
habitat selection in stochastic environments is largely lacking. Here, we take a step to 
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addressing this mathematical shortfall while at the same gaining new insights into the 
evolution of habitat selection for populations living in stochastic, patchy environments. 

the ideal-free distribution has been 



Since the classic paper Fretwell and Lucas 1969 



studied extensively from empirical, theoretical, and mathematical perspectives. Empir- 
ical support for ideal-free distributions exists for many taxa including fish |Godin and 



Keenleyside[ 


1984 


Oksanen et al. 


1995 Haugen et al. 


2006 


caster et al. 


1997 


, mammals [Beckmann and Berger 




2003 , 



For example, Oksanen et al. 1 1995 



birds | Harper 



1982 



Don- 



1995 



found that armored catfish in Panamanian stream 
pools were distributed such that the resource availability per catfish was equal in all 
occupied pools, despite significant variation in light availability across these occupied 
pools. Theoreticians have identified several "non-ideal" mechanisms (e.g. sedentarism, 
adaptive movement with finite speed, density-dependent dispersal) that, under equi- 
librium conditions, generate an ideal-free distribution Hastings, 1983, Cosner, 2005 



Gejji et al. , 2012 . For example, at equilibrium, sedentary populations achieve an ideal- 
free distribution provided, paradoxically, the populations initially occupied all habitat 
patches. While many early studies asserted that the ideal free distribution is an evo- 



lutionarily stable strategy (ESS) [Fretwell and Lucas 1969, van Baalen and Sabelis 



1993, Schreiber et al. 20001, only recent advanced nonlinear analyses fully verified this 



assert ion |Cressman et al. 2004, Cressman and Kfivan, 2006, 2010 Cantrell et al. 2007 



2010 20121. 



In nature, observed habitat occupancies are frequently less extreme than predicted 
by the ideal-free distribution: individuals underuse higher quality habitats and overuse 
lower quality habitats compared to theoretical predictions Milinski 1979, Tregenza 



1995 . Notably, populations occupying sink habitats have been documented in many 



species 



Sokurenko et al., 2006, Tittler et al., 2006, Robinson et al., 2008, Anderson 



and Geber 2010 . One possible explanation for these observations is that populations 
experience temporal as well as spatial variation in environmental conditions and, conse- 
quently, theory based on equilibrium assumptions tells an incomplete story. In support 
of this explanation, several theoretical studies have shown that occupation of sink habi- 
tats should evolve when temporal variation is sufficiently great in other habitats | Holt[ 



1997] |Jansen and Yoshimura[ [19981 l Holt and Barfield[ |2~00lj |Schreiberj [2012] . These 



theoretical developments, however, rely on linearizations of density-dependent models, 
and do not analyze the dynamics of competing genotypes, the ultimate basis for evo- 
lutionary change due to natural selection. Hence, these studies leave unanswered the 
question, "Does the linear analysis correctly identify competitive exclusion in pairwise 
interactions that is the basis for the analysis of evolutionarily stable strategies?" 

Within populations, individuals can exhibit different habitat selection strategies, and 



there is some evidence these differences can be genetically based [Via, 1990, Jaenike and 



Holt, 1991 . For instance, some individuals of the fruit fly species Drosophila tripunc- 



tata prefer tomato host plants (one potential habitat for its larvae) while others prefer 
mushrooms (another potential habitat), and these differences are based on two geneti- 



cally independent traits, settling behavior and ovipositor site preference Jaenike, 1985 



Jaenike and Holt |1991| found that genetic variation in habitat selection is common, espe- 



cially in arthropods and mollusks. Furthermore, they demonstrated using mathematical 
models that this genetic variation can stem from density-dependent regulation occurring 
locally within each habitat. Specifically, Jaenike and Holt write "frequency-dependent 
selection favors alleles that confer upon their carriers a preference for underused habitats, 
even if there is no genetic variation in how well individuals are adapted to the different 
habitat" [Jaenike and Holt 1991 p.S78l. Their analysis, however, doesn't account for 
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temporal fluctuations in environmental conditions and this raises the question, "Does 
environmental stochasticity facilitate or hinder the maintenance of genetic variation in 
habitat selection?" 

To answer the aforementioned questions, we provide an in-depth analysis of a model 
introduced in |Schreiber 2012 . The single genotype (i.e. monomorphic) version of this 



model and a characterization of its dynamics are given in Section 2. The competing 
genotype (i.e. dimorphic) version of the model and invasion rates of each genotype 
when rare are introduced in Section 3. In Section 4, we prove that these invasion rates 
determine the long-term fate of each of the genotypes. Specifically, if both genotypes 
have positive invasion rates when rare, then there is a positive stationary distribution 
under which the genotypes coexist. Alternatively, if one genotype has a negative inva- 
sion rate when rare, then it is asymptotically displaced by the other genotype. These 
result allows us to use the invasion rates when rare to explore conditions supporting a 
protected polymorphism for habitat selection. In Section 5, we refine the characteri- 
zation of an evolutionary stable strategy for habitat selection from iSchreiber 2012 in 



a mathematically rigorous manner, and provide an explicit formula for this ESS when 
there are two habitat types. Section 6 concludes with a discussion of how our results 
relate to the existing literature and identifies future challenges for the theory of habitat 
selection in stochastic environments. 

2. The Monomorphic Model 

To set the stage for two competing populations spread over several patches, we start 
with a single population living in one patch. Let Zt be the population abundance at 
time t > 0. The stochastic process (Z t )t>o is governed by the Ito stochastic logistic 
equation 

(2.1) dZt = Z t (n-KZt)dt + aZ t dW t , 

where /i is the intrinsic rate of growth of the population in the absence of stochasticity, n 
is the strength of intraspecific competition, a 2 > 0 is the infinitesimal variance parameter 
of the stochastic growth rate, and (Wt)t>o is a standard Brownian motion. The process 
(Zt)t>o is a strong Markov process with continuous paths. We call an object with such 
properties a diffusion. 

As shown in our first proposition, the process {Zt)t>o lives in the positive half line 
M++ := (0, oo); that is, if we start it in a strictly positive state, then it never hits zero. 
Furthermore, the long-term behavior of the process is determined by the stochastic rate 

2 

of growth \i— \. When the stochastic growth rate is negative the population abundance 
converges asymptotically to zero with probability one. On the other hand, when this 
parameter is positive the distribution of the abundance converges to an equilibrium 
given by a Gamma distribution. These results are well-known, but, as introduction to 
the methods used to prove our main results, we provide a proof in Appendix A. 

Proposition 2.1. Consider the diffusion process (Zt)t>o given by the stochastic differ- 



ential equation (2.1). 

• The stochastic differential equation has a unique strong solution that is defined 
for all t > 0 and is given by 

z = Z 0 exp(( / u-a 2 /2)t + t jW f ) 

1 + Z 0 ^ Jj exp(Gu - o 2 /2)s + aW s )ds 

• // Zq = z > 0, then Zt > 0 for all t > 0 almost surely. 
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// fi — %2 < 0, then limi_s.oo Zt = 0 almost surely. 

V f 1 ~ %~ = 0; then liminf^oo Z t = 0 almost surely, limsup^^ Z t = oo almost 
surely, and lim^oo j J 0 * Z s ds = 0 almost surely. 

2 

If fj, — > 0; f/ien (Zt)t>o has a unique stationary distribution p on M++ wif/i 
Gamma density g(x) = r ^ ek x k ~ 1 e~» , where 



— and k '.= — „ — 1. 
2k ct 2 



Moreover, if Zq = z > 0, i/ien 



1 



o 



lim - / h(Z s ) ds = h(x)g(x) dx almost surely 

t—>-oo t 



for any B or el function h : R++ — > K wii/i / 0 °° < oo. In particular, 

If* If a 2 \ 
lim -/ Z s ds = — ■ [a almost surely. 

t^-co t J 0 k \ 2 J 

Next, we consider a population living in a spatially heterogeneous environment with 
n different patches. These patches may represent distinct habitats, patches of the same 
habitat type, or combinations thereof. The abundance of the population in the i-th 
patch at time t > 0 is X\. Let (Xl)t>o be given by 

(2.2) dX\ = X\ (m - KiX'l) dt + X\ dE l t , 

where \x% is the intrinsic rate of growth the population in patch i in the absence of stochas- 
ticity, Ki is the strength of intraspecific competition in patch i, and El = Yljlji^t f° r 
a standard multivariate Brownian motion (B 1 , . . . ,B n ) T on W 1 and an n x n matrix 
r := (7ij). The infinitesimal covariance matrix for the non-standard Brownian motion 
(£/,..., is S = (ay) := T T T. 



The populations in the various patches described by equation (2.2) are coupled only by 
the spatial correlations present in the driving Brownian motion (E£, . . . , Ef). We further 
couple the population dynamics across patches by assuming the fraction of population 
in patch % equals a.% for all time. This spatial distribution can be realized at the scale of 



the individual when, as described in greater detail in Remark 2.2, individuals disperse 
rapidly and independently of one another in such a manner that the fraction of time 
spent in patch i equals oti for each individual. Under this assumption, we call a = 
(ai,Q!2, • • • ,ct n ), with a.i > 0 for all 1 < i < n and Y17=i oi% = 1, a patch-selection 
strategy. Continuing to denote the abundance of the population in the i-th patch at 
time t > 0 as XI, we have X\ = oiiXt, where Xt = Y17=i^-t ^ s the total population 
abundance at time t > 0. If we impose these constraints on (X 1 ...,X n ), then it 
is heuristically reasonable that the process X is an autonomous Markov process that 
satisfies the SDE 

n n 

(2.3) dX t = X t ^2ai (m - KiaiXt) dt + X^ (*i dE\. 



Remark 2.2. One way to justify the formulation of ( 2.3 ) rigorously is to first modify ( 2.2 ) 
to obtain a system of SDEs explicitly accounting for dispersal. Suppose that individuals 
disperse from patch i to patch j at a rate 5dij for some fixed rate matrix D = (dij). As 
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usual, we adopt the convention da = — d ij- The resulting system of SDEs is 

(2.4) dX\ = X\ (jn - KiXf) dt + S^ Xfdji dt + X\dE[. 

j 

Assume that the rate matrix D has a unique stationary distribution a; that is, ay > 0 
for l<j<n, Y!i=i a i = 

n 

(2.5) J>^ = 0 

i=i 

for 1 < i < n. In this case, a vector (y 1 , . . . , y n ) satisfies 



(2.6) 

for 1 < % < n if and only if 
(2.7) 



y J = ocjC 



for 1 < j < n for some constant c. Moreover, summing (2.7) we find that 

n 



(2. 



Note that by (2.5) we can write the drift term in (|2.4|) that contains 5 as 
(2.9) 



5 ^ X t d ji dt = 6 J2^t - a j X t) d ji dt 



where X t := Y% =1 X\. using §T?§ and we see that (as 1 , . . . , x n ) and x := £? =1 x i 

are such that 



(2.10) 

for i = 1, . . . , n if and only if 
(2.11) 



^^(x- 7 — ctjx)dji = 0 



5> 



1=1 



for 1 < j <n. 



It follows from (2.9) and the equivalence between (2.10) and (2.11) that as 5 increases 



the solution of (2.4) experiences an increasingly strong drift towards the one-dimensional 
subspace 

{(xi, . . . ,x n ) : Xi = cti(xi H h x n ), i = 1, . . . , n}. 



In the limit S — > oo, it is plausible that the system (2.4) converges to one for which 



X\ = atiX t , 

where Xt := X\ + • • • + X™, and the total population size Xt satisfies the autonomous 



one-dimensionl SDE (2.3) with Xt replaced by Xt. This heuristic for obtaining (2.3) as 



a high dispersal rate limit of (2.4) can be made rigorous by applying Theorem 6.1 from 
Katzenberger |1991 . 
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Let x • y = Y17=i x iVi denote the standard Euclidean inner product and define another 
inner product (•, -) K by (x,y) K := J27=i K i x iVi- Since (a ■ E t )t>o is a Brownian motion 
with infinitesimal variance parameter a • Ea, ( |2.3[ ) can be expressed more simply as 

(2.12) dX t = X t (a ■ fx - (a, a) K X t ) dt + X t Va ■ ZadW t , 

where Wj is a standard Brownian motion. 



The total population (Xt)t>o defined by (2.12) behaves exactly like the one-patch 
case defined by ( |2.1[ ) with the parameters fx — >• fx ■ a, k — > (a,a) K and a — > V a ■ E<x 
In particular, {Xt)t>o is a diffusion process and we have the following immediate conse- 



quence of Proposition 2.1 



Proposition 2.3. Consider the diffusion process (X)t>o given by (2.12). 

• // Xq = x > 0, then Xt > 0 for all t > 0 almost surely. 

• If a ■ [i — < 0, then lim^oo X t = 0 almost surely. 

• If a ■ fi — ot '^ ot = 0, then liminft^oo Xt = 0 almost surely, limsup t _^ 00 Xt = oo 

l ft 



almost surely, and lim^oo j J Q X s ds = 0 almost surely. 

2 

r(fc)6"= 



If a- fi — 9L r^ L > 0, then the process (Xt)t>o has a unique stationary distribution 
Px on M++ with Gamma density g(x) = r ,A p x fc ~ 1 e~g , where 



a ■ Ea 2a- fx 

and k := — 1. 



2(a, a) K a • Ea 

Moreover, 



t /-oo 



lim - / h(X s ) ds = I h(x)g(x) dx almost surely 

t^co t J 0 Jq 

for any B or el function h : M++ — > M with \h(x)\g(x)dx < oo. In particular, 
lim - / X s ds = a ■ ( u | almost surely. 



For the dynamics ( 2.2 ) in patch i, Proposition 2.1 implies that if there was no coupling 
between patches by dispersal, then then population abundance in patch i would converge 
to 0 if fii — a a 1 2 < 0 and converge to a non-trivial equilibrium if fn — cru/2 > 0. As 



noted by Schreiber |2012 and illustrated below, the spatially coupled model is such that 



the population can persist and converge to an equilibrium even when m — a a/ 2 < 0 for 
all patches. 

Persistence of coupled sink populations in symmetric landscapes. Consider a 
highly symmetric landscape where fi% = r, an = a 2 > 0 for all i , Ki = a for all i, 
and Oij = 0 for all i ^ j. If individuals are equally distributed across the landscape 
(oj = 1/n for all i), then 

an a 2 a ■ Ea a 2 

fj,i = r and a ■ fx 



2 2 r 2 2n 

The increase in the stochastic growth rate from r — a 2 /2 for an isolated population to 
r — a 2 j (2n) for the spatially coupled population stems from individuals spending equal 
time in patches with uncorrelated environmental fluctuations. Specifically, the environ- 
mental variance experienced by individuals distributing their time equally amongst n 
uncorrelated patches is n times smaller than the environmental variance experienced by 
an individual spending their time entirely in one patch. Whenever a 2 > 2r > a 2 /n, this 
reduction in variance allows the entire population to persist despite patches, in and of 
themselves, not supporting population growth. 
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3. Dimorphic model and invasion rates 

To understand the evolution of patch-selection strategies, we now consider competi- 
tion between populations that only differ in their patch-selection strategy. Let Xt and 
Yt be the total population sizes at time t > 0 of two populations playing the respective 
patch selection strategies a = (ati, a-i, • • • , ot n ) and (3 = (Pi,P2, ■ . ■ , P n ), so that the den- 
sities of the populations in patch i are ctiXt and PiYt at time t > 0. The dynamics of 
these two strategies are described by the pair of stochastic differential equations 



(3.1) 



Since 



dX t = Xt^ai (m - Ki(ociX t + PiY t )) dt + X t ^2 cm, dE. 

i=l i=\ 
n n 

dY t = Y t ^2f3i {fjn - Ki(aiX t + PiY t )) dt + Y^ Pi dE\. 



d[X,X] t = X?a-T,adt 
d[YY] t = Y t 2 (3-Z(3dt 
d[X,Y] t = X t Y t a-Zpdt, 

the diffusion process ((X t ,Y t )) t >Q for the spatially coupled, competing strategies can be 
represented more compactly as 

dX t = X t [fj, ■ a — (a, (3) K Y t — (a, a) K X t ] dt + X t \l a ■ Sa dUt 

(3 ' 2) dY t = Yt[fi-p- (a, P) K X t - (P, P) K Y t ] dt + Yty/P ■ £/? dV u 

where (U, V) is a (non-standard) Brownian motion with covariance structure d[U, U]t = 
dt, d\V,V]t = dt, and d\U,V\t = ^^ S/ L= = dt. Using a construction similar from 



Remark 2.2, system (I3.2I) can be seen as a high dispersal limit. This system exhibits 



a degeneracy when U = V i.e. , "' S ^L= ; = 1. If S is nonsingular, then, by the 
Cauchy-Schwarz inequality, this degeneracy only occurs if a = p. We do not consider 
this possibility in what follows. 

To determine whether the two populations coexist or one displaces the other, we 
introduce the invasion rate Z(a, P) of a population playing strategy P when introduced 
at small densities into a resident population playing strategy a. As shown in the next 
proposition, this invasion rate is defined by linearizing the dynamics of Y and computing 
the long-term population growth rate I(a, P) associated with this linearization. When 
X(a, P) > 0, the population playing strategy P tends to increase when rare. When 
I(a, P) < 0, the population playing strategy p tends to decrease when rare. 

Proposition 3.1. Consider the partially linearized system 

dXt = Xt\u-a — (a, a) K Xt\ dt + XtV a ■ Sa dUt 

(3.3) A A L _ A , 

dY t = Y t [n-P- (a, P) K X t ] dt + Yty/p-XpdVt. 

Assume Xq > 0 and Yq > 0. 

If a ■ (jJt — Sa/2) > 0, so the Markov process X has a stationary distribution concen- 
trated on then the limit lim^oo lo ^ y * exists almost surely and is given by 

(3.4) X(a, P) = P ■ ( M - S/3/2) - ■ (jjl - Ea/2) . 

(a, a) K 
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On the other hand, if a-([i — Sa/2) < 0, so that lim^oo j X s ds = 0 almost surely, 
then the limit lim^oo l ° e t Yt exists almost surely and is given by 

(3.5) X(a,/3) = /3-( / u-S/3/2). 

Proof. By Ito's lemma, 

dlogY t = ((j,-B-{a,B) K X t ) dt + y/p-XpdVt-^(P-XP)dt. 



2.3 



Assume that \l- a — - 2 - > 0. By Proposition 

If 1 - a 
lim - / X s ds = — • (2fi — Sa) almost surely. 

t J Q 2(a,a) K 



Therefore, 



lim — - = 8 ■ (j, ^—^ a ■ (2u — Sa) /3 ■ S/3 almost surely, 

t-5-oo t 2{a,a) K 2 



as claimed. 



On the other hand, assume that /j, • a — < 0. By Proposition 



2.3 



1 f l - 

lim - / X s ds = 0 almost surely. 

t— >oo t 



0 



Therefore, 



lim — - — * = a ■ (3 B ■ TiB almost surely, 

t— >oo t 2 



again as claimed. □ 

In the next proposition, we show that if a population playing strategy B cannot invade 
a population playing strategy a (i.e. I(a,B) < 0), then the population strategy a can 
invade the population playing strategy 8 (i.e. I(B,a) > 0). This suggests, as we will 
show in the next section, that such a strategy a should exclude strategy B. 

Proposition 3.2. Suppose that a- (fi — Sa/ '2) > 0 and Z(a, 8) < 0. Then, 1(8, a) > 0. 

Proof Set A := a ■ (/Jb- Ea/2) and B := B ■ (fj, - T.8/2). Assume that A > 0 and 
I (a, 8) < 0. To show that 1(8, a) > 0, we consider two cases, B < 0 and B > 0. 
Suppose B < 0. Then, 1(8, a) = A > 0 by Proposition |3. 1| and by assumption. 
Alternatively, suppose that B > 0. Then 

l(a, B) = B- A { ^^ and 1(8, a) = A - B { ^^ 
{a,a) K {P,B) K 



by Proposition 3.1 Assume, contrary to our claim, that Z(a,8) < 0 and T(B,a) < 0. 

1/2 1/2 

From the Cauchy-Schwarz inequality (x,y) K < (x,x) K (y,y)n we get 

B( a ,a) l J 2 (B,B) l J 2 > B(a, B) K > A(8,8) K 
^(a,a)y 2 (/3,/3)y 2 > A(a,8) K > B(a,a) K . 

1/2 1/2 1/2 

The above inequalities yield the contradiction B(a, a) K >A(B,B) K and B(a,a) K < 
A{8,8) l J 2 . □ 
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An immediate consequence of Proposition 3.1 is the following corollary. This corollary 



implies that if a population playing strategy f3 can invade a population playing strategy 
a and a population playing strategy a can invade a population playing strategy /3, then 
a single population playing strategy a converges to a non-trivial equilibrium and the 
same is true of a single population playing strategy (5. This suggests, as we show in the 
next section, that under these conditions these two strategies should coexist. 

Corollary 3.3. The invasion rate satisfies I(a,f3) < f3 ■ (fj, — E/3/2). In particular, if 
T(a, (3)>0 and J(/3, a) > 0, then a • (// - Ea/2) > 0 and (3 ■ (fx - E/3/2) > 0. 

4. Exclusion and protected polymorphisms 

Our main results about the dimorphic process (X, Y) is that the invasion rates de- 
termine the long-term fate of competing strategies. If the invasion rates predict that 
strategy j3 cannot invade a population playing strategy a, then the population play- 
ing strategy a drives the population playing strategy (3 asymptotically to extinction as 
shown in the following theorem. We give a proof in Appendix B. 

Theorem 4.1. If a • (p - Ea/2) > 0 and l(a,(3) < 0, then, for (x,y) G the 
probability measures 

rt 



\ f p^y){( Xs ,Y s )e-}ds 
t Jo 



converge weakly as t —> oo to p x (g) So, where p x is the unique stationary distribution of 
Xt concentrated on and 5$ is the point mass at 0. 

On the other hand, if the invasion rates predict that each strategy can invade when 
rare, then the following theorem proves that the competing strategies coexist: for any 
initial conditions the joint distribution of (Xt,Yt) converges as t — > oo to a probability 
distribution on with density ip and, moreover, for any Borel set B C K+ + the long 
term proportion of times t for which (Xt, Yt) spends in B converges to 



jjj(x,y) dxdy. 



B 



A proof is given in Appendix C. In order to appreciate the assumptions of the theorem, 
it helps to recall Corollary |3.3| which says that if T(a, f3) > 0 and X(/3, a) > 0 then 
a ■ (n — Ea/2) > 0 and (3 ■ (fi — S/3/2) > 0 so that a single population playing strategy 
a or f3 will persist. 

Theorem 4.2. Suppose that I (a, /3) > 0 and I((3,a) > 0. Then, there exists a unique 
stationary distribution it of (X, Y) on that is absolutely continuous with respect to 
Lebesgue measure. Moreover, for any bounded, measurable function f : — > M. 

(4.1) lim - / f(X s ,Y s )ds= [ f(x,y)ir(dx,dy) almost surely. 

t^co t J 0 J M 2 ++ 

Furthermore, the process (X, Y) is strongly ergodic, so that for any initial distribution 
q one has 

(4.2) lim d TV (F q {(X t , Y t ) e •}, vr) = 0, 

t— >oa 

where dxv is the total variation distance. 
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From the perspective of population genetics, the coexistence of these two strategies 
corresponds to a protected polymorphism: each strategy (a morph) increases when 
rare and, therefore, is protected from extinction. This protection from extinction, how- 
ever, is only ensured over ecological time scales as mutations may result in new morphs 
that can displace one or both coexisting morphs Ravigne et al. 2004| . The concept of 
protected polymorphisms was introduced by Prout 1968] when studying deterministic 



models of competing haploid populations in a spatially heterogenous with overlapping 
generations. Turelli et al. [2001 extended this concept to stochastic difference equa- 
tions for competing haploid populations with a constant population size. Theorem 4.2 
provides a mathematically rigorous characterization of protected polymorphisms for our 
stochastic models with fluctuating population sizes. 

Theorem |4.2| implies that coexistence depends on the intrinsic stochastic growth rate 
of the populations and the competitive effect of each population on the other. The 
intrinsic stochastic growth rates are given by 

r a = a ■ (fi — Sq/2) and rp = (3 • — S/3/2). 

While the competitive effect of the population with strategy a on the population with 
strategy f3 is given by the ratio of the magnitude of a projected in the ft direction 
(i.e. {p/\\p\\ K ,a) K where \\P\\ K = y/{P,P) K ) divided by the magnitude of p (i.e. \\p\\ K ). 
Mathematically, the competitive effect of a on f3 and the competitive effect of j3 on a 
are given by 

(P/\\P\\ K ,a) K WIHU/3). 
' " \\a L 



Provided r a and rp are positive, Theorem 4.2 implies that there is a protected polymor- 
phism if 



(4.3) 



Tg 

> C Bn and — > C, 



Q,/3- 



In words, the relative intrinsic stochastic growth rate of each population must exceed 
the competitive effect on itself due to the other population. Conversely, if one of the 



inequalities in (4.3) is reversed, then Theorem 4.1 implies that one population excludes 



the other. Unlike the standard Lotka-Volterra competition equations, Proposition 3.2 



implies that both inequalities in (4.3) cannot be simultaneously reversed and, conse- 



quently, bistable dynamics are impossible. 

Environmental stochasticity impedes protected polymorphisms in symmetric 
landscapes. Consider a landscape where all patches have the same carrying capacities 
(e.g. Ki = 1 for all i), the same intrinsic rates of growth (i.e. \i{ = a for all i), and the 
same amount of uncorrelated environmental stochasticity (e.g. an = a 2 for all i and 



<Tjj = 0 for i 7^ j). Then the protected polymorphism inequalities (4.3) become 



(4.4) 



a — a 



I«II72 



a — o i 



2 /2 



> C/3,q and 



a 



72 



\a\\ 2 /2 



>a 



/2 



where the only a dependency is on the left hand sides of both inequalities. As r _ a 2\\a\\y2 
is a decreasing function of a 2 whenever |^| < 1 and an increasing function of a 2 when- 



ever ^ 
phism 



11/ 

> 1, it follows that the set of set of strategies supporting a protected polymor- 
A(a 2 ) = {(a,P) : (144} holds} 
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Figure 1. Protected polymorphisms and exclusion in two-patch land- 
scapes. Contour plots of I(a, (3) where lighter shades correspond to 
higher values of I(a,f3). The regions where I(a, /3)l(f3,a) > 0 are de- 
lineated by the solid curves and correspond to parameter combinations 
supporting a protected polymorphism. Regions where I(a, /3)I(a, (3) < 0 
correspond to strategies that cannot coexist. The dashed-dotted and dot- 
ted curves indicate how regions of coexistence and exclusion change for 
higher and lower levels of environmental stochasticity a 2 , respectively. 
In panel A, the landscape is spatially homogeneous with \i = (1, 1), 
k = (1, 1) and E = a 2 I where / is the 2x2 identity matrix. In B, 
the landscape is spatial heterogeneous with respect to the deterministic 
carrying capacities k = (3, 1) and the remaining parameters as A. 



is a decreasing function of a 2 i.e A((J2) is a proper subset of A(a 2 ) whenever 02 > <J\ > 0. 
Figure [l]A illustrates this conclusion for a two-patch landscape. Intuitively, increas- 
ing environmental stochasticity in these symmetric landscapes reduces the stochastic 
growth rate for all strategies and, thereby, makes it less likely for populations to persist 
let alone coexist. For asymmetric landscapes, how the set A(a 2 ) of protected poly- 
morphisms varies with a 2 is more subtle, as illustrated in Figure [lp. In this case, some 
protected polymorphisms are facilitated by environmental stochasticity, while other pro- 
tected polymorphisms are disrupted by environmental stochasticity. 

For the symmetric landscapes, we can identify a strategy that displaces all others. 
Namely, the strategy a = . . . , -) of visiting all patches with equal frequency. This 
strategy maximizes the function function a 1— > a — a 2 \\a\\ 2 /2. Hence, if we consider a 
competing strategy f3 7^ a, then a-/3=||a|| 2 = ^ and 



Z{a^) = a-a 2 W\\ 2 /2 
= a-a 2 \\P\\ 2 /2 



a- 15 



a 



a 2 \\a\\ 2 



a 2 \\a\\ 2 



< 0. 
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e.g. the invasion rates are negative along the vertical transect ot\ = 1/2 in Figure. [T]A.. 
This strategy a is an example of an evolutionarily stable strategy that we discuss further 
in the next section. 



5. Evolutionarily stable strategies 



The concept of an evolutionary stable strategy was introduced by Maynard Smith 



and Price 1973| . Loosely stated, an evolutionary strategy is a strategy that cannot 



be invaded by any other strategy and, consequently, can be viewed as an evolutionary 
endpoint. For our models, we say patch selection strategy a is an evolutionarily stable 
strategy (ESS) if X(a,/3) < 0 for all strategies f3 ^ a. In light of Theorem 4.1, an 
ESS not only resists invasion attempts by all other strategies, but can displace all other 
strategies. An ESS a is called a pure ESS if cti = 1 for some patch i, otherwise it 
is a mixed ESS. Our next result provides an algebraic characterization of mixed and 
pure ESSs. However, it remains to be understood whether these ESSs can be reached 



by small mutational steps in the strategy space (i.e. are convergently stable Geritz 
etaL|[l997| ). 



Theorem 5.1. Assume that the covariance matrix S is positive definite and that there 
is at least one patch selection strategy which persists in the absence of competition with 
another strategy; that is, that max a a ■ (fj, — T,a/2) > 0. 

Mixed strategy: An ESS a with ai > 0 for i £ I with I C {1,2,. . . , n} and 
\I\ > 2 satisfies 



a • Ea a ■ (2u - Sa) 

5-1 0 — =m-Kiai — — > 

2 2(a,a) K -f-' 



(TijCtj 



for alii E /. Conversely, if\I\ = n, then a strategy a satisfying (5.1) is an ESS. 



Pure strategy: The strategy on = 1 and aj = 0 for j ^ i is an ESS if and only if 



(5-2) N- a f<- a f + °H- a i 



for all j i. 

Furthermore, in the case of n = 2, there exists a mixed ESS whenever the reversed 
inequalities 

(5-3) /', ? ? • a " 



2 J 2 



hold for i = 1, 2 and j ^ i. 



The first statement of Theorem |5.1| provides a sufficient and necessary condition for 
a mixed ESS utilizing all patches. For example, in a symmetric landscape (as described 
in the previous section), this ESS condition is only satisfied for a = (1/n, 1/n, . . . , 1/n). 



The second statement of Theorem 5.1 provides a characterization of when using only 



a single patch is an ESS. Since the right hand side of equation (5.2) is negative, using 
patch i can only be an ESS if all other patches have a negative stochastic rate of growth, 
lij — < 0 for all j i. However, even if only patch i has a positive stochastic growth 
rate, an ESS may use the other patches, as we illustrate next for two-patch landscapes. 
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Figure 2. ESS for patch selection (left) and mean population abun- 
dance (right) in a source-sink landscape. Parameter values: n = 2, 



0-22 = (Ji2 = 0, /i = (l,/x 2 ), and At = (1, 1). 



ESSs in two-patch, uncorrelated landscapes. For an uncorrelated two patches 
landscape (i.e. n 
whenever 



2 and a\i = 0), Theorem 5.1 implies that there is a mixed ESS 
Ml > — CT22/2 and \i2 > —an/2 



Hi + a ■ Sa 



(5.4) 

and this ESS satisfies 

(5.5) Qj 

Ki{H ■ a - a - Ha/2) /{a, a) K + an 

Equation ( |5.4[ ) implies that even if deterministic growth in patch 2 is strictly negative 
(i.e. H2 < 0), then there is selection for movement into this patch provided the variance 
of the fluctuations in patch 1 are sufficiently large relative to the intrinsic rate of decline 
in patch 2 (Fig. [2}. 

In the limit of no noise (i.e. an 4- 0 for i = 1, 2), equation (5.5) becomes 

a- - ^ 

Kifi ■ a/(a,a) K 

While our results do not apply to the deterministic case, this limiting expression for the 
ESS suggests, correctly, that the ESS for the deterministic model satisfies 

Hi/ K i 



whenever /ij > 0. 



In other words, the fraction of individuals selecting patch i is proportional to the equi- 



librium density Hi/ K i supported by patch i. Equation (5.5) implies that adding stochas- 
ticity in equal amounts to all patches (i.e. an = a 2 for all i) results in an ESS where, 
relative to the deterministic ESS, fewer individuals select patches supporting the high- 
est mean population abundance and more individuals selecting patches supporting lower 
mean population abundances (Fig. [3]). 



6. Discussion 



Habitat selection by organisms is a complex process determined by a mixture of 
genetic, developmental, ecological, and environmental factors. For ecologists, habitat 
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Figure 3. The effect of the deterministic carrying capacities and envi- 
ronmental stochasticity on the ESS for patch selection (left) and mean 
population abundance (right) in a two-patch landscape. The ratio 
ftl/(«2 + K 2) corresponds to the ratio of the deterministic carrying ca- 
pacity (^2/^2) m patch 2 to the sum of the deterministic carrying ca- 
pacities (/U1/K2 + ^2/1^2) when fi\ = ^2 = 1- Parameter values: n = 2, 
en = C22 = c 2 , C12 = 0, n = (1, 1), and « = (1, « 2 )- 



selection plays a fundamental role in determining the spatial and temporal distribution 
of a population Rosenzweig, 1981, Orians and Wittenberger 1991] . For evolutionary 



biologists, habitat selection determines the suite of environmental factors driving local 
adaptation [Edelaar and Bolnick 2012 . Indeed, in the words of the eminent evolution- 
ary biologist Ernst Mayr, "With habitat and food selection - behavioral phenomena - 
playing a major role in the shift into new adaptive zones, the importance of behavior 

[Mayr[ |1963[ p. 604]. 



Here, we 



in initiating new evolutionary events is self-evident" 
examined how spatial and temporal heterogeneity in demographic rates across multiple 
habitat patches influence the dynamics of competing populations who only differ in their 
habitat patch selection preferences. We assume that habitat selection has a genetic ba- 
sis (e.g. genes that influence the physiological or neurological capacity of individuals to 
detect and respond to habitat cues) and that genetic differences in habitat choice have 
no pleiotropic effects on habitat specific fitness. Our analysis reveals that, generically, 
only two outcomes are possible, coexistence or displacement of one population by the 
other for all initial conditions, and that these outcomes are determined by the invasion 
rates of populations when rare. In addition to providing a mathematically rigorous jus- 
tification of prior work, our analysis provides new insights into protected polymorphisms 
for habitat selection and raises several questions about evolutionary stable strategies for 
habitat selection. 

Protected polymorphisms correspond to populations of competing genotypes exhibit- 



ing negative frequency-dependence: each population tends to increase when rare Prout 



1968 . In the case of patch selection, these competing populations differ in the frequen- 
cies in which they select habitat patches. In a survey of the empirical literature, Jaenike| 
and Holt 1991 found "that genetic variation for habitat selection is common, especially 
in arthropods and mollusks, the groups that have been studied most frequently." More- 
over, they argued that some of this variation may be maintained through protected 
polymorphism. Specifically, "in a haploid model without intrinsic fitness differences 



Downloaded from http://biorxiv.org/on September 18, 2014 



EVOLUTION OF PATCH-SELECTION IN STOCHASTIC ENVIRONMENTS 15 

among genotypes [i.e. soft selection], genetic variation in fixed habitat preferences may 



be maintained stably" [Jaenike and Holt 1991, pg. S83]. We provide a general ana- 



lytic criterion (see, inequality (4.3)) characterizing these protected polymorphisms for 
spatially and temporally variable environments. This criterion depends on the intrinsic 
fitnesses (r Q and r«) of each population and their competitive coefficients (C a p and 
Cp a ) that characterize the effect of each population on the other. Competitive effects 
are greatest when there is an overlap in patch use and one population tends to select the 
patches with the higher carrying capacities more than the other population. Intuitively, 
by occupying patches with a larger carrying capacities, populations achieve higher re- 
gional densities. Coupled with overlap in patch use, these higher densities result in a 
greater competitive impact of one population on another. A protected polymorphism 
occurs when the relative fitness of each population (e.g. r a /rp for strategy a) is greater 
than the competitive effect of the other population on it (e.g. r a /rp > C^ a for the 



population playing strategy a). Hence, as in the case of species coexistence |Chesson 



2000 , protected polymorphism are most likely when fitness differences are small (i.e. 
r a /rp « 1) and competitive effects are small (i.e. both C a ^ and Cp jCl < 1). Envi- 
ronmental stochasticity solely effects the intrinsic fitness terms and can facilitate or 
inhibit protected polymorphisms. For landscapes in which all patches experience the 
same degree of uncorrelated, temporal variation, environmental stochasticity has an 
inhibitory effect as it magnifies fitness differences between competing strategies (e.g. 
r a/ r /3 increases with environmental stochasticity). For asymmetric landscapes, how- 
ever, temporal variability can facilitate polymorphisms by reducing fitness differences of 
competing strategies. 

In contrast to protected polymorphisms, our analysis reveals that populations play- 
ing an evolutionarily stable strategy (ESS) for patch selection not only thwart invasion 
attempts by all other strategies but also can invade and displace a population playing 
any other strategy. Furthermore, our analysis provides a mathematically rigorous jus- 



tification of an earlier characterization of ESSs Schreiber, 2012 . This characterization 



implies that populations playing the ESS always occupy source habitats (i.e. patches 
where \xi — ofj/2 > 0). Indeed, consider a population playing strategy a that does not 
occupy some source patch, say patch i. Then a different behavioral genotype (3 that 
only selects patch i can invade as I(a,/3) = fii — af/2 > 0. In the limiting case of a 
deterministic environment, our characterization of the ESS recovers the classic result of 



McPeek and Holt 1992 : the fraction of time spent in a patch is proportional to the 
carrying capacity of the patch. Adding environmental stochasticity generally results in 
populations playing the ESS decreasing the time spent in the patches with larger carry- 
ing capacities and possibly making use of sink patches (i.e. patches where Hi—af/2 < 0). 
This shift in patch choice can be viewed as a spatial form of bet hedging: individuals 
increase fitness by decreasing the variance in their stochastic growth rate at the expense 



of their mean growth rate Childs et al. ,2010 



We are able to show that for two patch landscapes there always exists an ESS for 
patch selection. However, several questions remain unanswered. First, what happens 
for landscapes with more than two patches? Is there always an ESS? Second, while we 
know that a population playing an ESS can displace a monomoprhic population playing 
a different strategy, can it displace polymorphic populations? Finally, are ESSs always 



convergently stable Geritz et al. 19971? If there are positive answers to this final suite 



of questions, then ESSs can be generally viewed as the ultimate evolutionary end state 
for patch selection strategies. 
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Going beyond the models considered here, studying the evolution habitat use faces 
many challenges. Our models assume that populations spend a fixed fraction of time 
in each patch and do so instantaneously. What happens if we relax these assumptions? 
For example, if populations are more ideal and able to track changes in population 
density instantaneously, then we have something closer to the classical notion of ideal 
free movement Fretwell and Lucas, 1969 . For these populations, what is the optimal 
(in an evolutionary sense) density-dependent strategy? Moreover, can such a strategy 
displace the static strategies considered here? Alternatively, if populations are less ideal 
and diffusing randomly on the landscape, what happens then? The linear version of 
this question was tackled in part by Evans et al. 2013 . However, the mathematical 



analysis for analogous stochastic models with density-dependent feedbacks is largely 
unexplored. Going beyond single species, the coevolution of patch selection among 
interacting species has a rich history for spatially heterogeneous, but temporally ho- 



mogeneous environments van Baalen and Sabelis 


1993, 


Kf ivan 


1997 


Schreiber et al. 


2000, van Baalen et al. 2001, Schreiber et al. 2002, Cressman et al. 


2004, Schreiber 


and Vejdani 




2006 




Cantrell et al. 




2007 . For example, spatial heterogeneity can se- 



lect for the evolution of contrary choices in which the prey prefers low quality patches 
to escape the predator and the predator prefers high quality patches to capture higher 



quality food items (Fox and Eisenbach, 1992, Schreiber et al. 2000 . Understanding how 



environmental stochastic influences this coevolution of patch choice and the community 
level consequences of these coevolutionary outcomes provides a plethora of important, 
yet largely untouched challenges for future work. 



Appendix A: Proof of Proposition 12.11 
The stochastic differential equation for Z is of the form 
(6.1) dZ t = b(Z t ) dt + a(Z t ) dW t , 

where b(z) := fiz — kz 2 and cr(z) := az. It follows from Ito's existence and uniqueness 
theorem for strong solutions of stochastic differential equations that this equation has a 
unique strong solution up to possibly a finite but strictly positive explosion time. 
Set Rt := log Zt for t > 0. By Ito's lemma, 



3.2) 



dRt 



ft 



(7 

~2 



k exp(R t ) ] dt + a dW t . 



It follows from the comparison principle of Ikeda and Watanabe (see Chapter VI The- 
orem 1.1 of Ikeda and Watanabe [1989] ) , Theorem 1.4 of Le Gall |1983|, or Theorem 
V.43.1 of |Rogers and Williams| [2000] ) that 



(6.3) 



Rt<Ra+[v 



<j- 



t + aW t . 



and so Z does not explode to +oo in finite time. Moreover, since r 4 |J - ne r is a 
bounded, uniformly Lipschitz function on (— oo,0] it follows from Ito's existence and 
uniqueness theorem that R does not explode to — oo in finite time, so that Z does 
not hit 0 in finite time. We could have also established this result by using the scale 
function and speed measure calculated below to check Feller's necessary and sufficient 
for the boundary point of a one-dimensional diffusion to be inaccessible - see Theorem 
23.12 of |Kallenberg] (2002 . 
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IT 



It is not hard to check using Ito's lemma that an explicit solution of the SDE is 



Z t . 



Z 0 exp((n-a 2 /2)t + aWt 



l + Z 0 ^ £ exp((/x - a 2 /2)s + aW s ) ds 
We see from the inequality (6.3) that if jj, — a 2 /2 < 0, then lim, 



t^oo Z t = 0 almost 
surely. 

We use the theory based on the scale function and speed measure of a one-dimensional 



diffusion (see, for example, Chapter 23 of Kallenberg 2002| or Sections V.6-7 of Rogers 
and Williams |2000| ) below to establish that Z is positive recurrent with a unique sta- 
tionary distribution when fi—a 2 /2 > 0. Similar calculations show that Z is null recurrent 
when n — a 2 /2 = 0, and hence liminft_ s . 00 Zt = 0 almost surely and limsup^oo Zt = oo. 



It follows from (6.2) and the comparison principle that if Z' and Z" are two solutions 



of (6.1) with respective parameters fj,',K',cr' and //',«/', <t" satisfying // < /i" , 



a' = a" and the same initial conditions, then Z[ < Z" . We will show below that 



lim - 

t— >oo t 



Ze ds 



1 

K 



(v-a 2 /2) 



almost surely when \i — a 2 / 2 > 0, and hence 



lim - 

t— >oo t 



Ze ds = 0 



almost surely when fi — a 2 /2 = 0. 

We now identify the scale function and speed measure of the one-dimensional diffusion 
Z. A choice for the scale function is 

t 2b{z) 



six 



(6.4) 



cxp 



o*(z) 



dz dy 



2K {y-a) 



dy 



for arbitrary a, c G (recall that the scale function is only defined up to affine 

transformations). If we set a = (as') o then 

ds(Zt) = a(s(Z t )) dW t 

and the diffusion process s(Z) is in natural scale on the state space s(M ++ ) with speed 
measure m that has density \. 

The total mass of the speed measure is 

1 . f 1 



m(R ++ ) 



_) ° 2 (x) 



dx 



1 



(6.5) 



a 2 a 



{auf 



2^/(7 



-2ti/a 2 



H{u-a) 



) {(<rJ) ° a- 
du 



i) 2 (x; 



dx 



1 



(u)s'(u) 



du 



\ a , 
oo 



^-2 



du. 



By Theorem 23.15 of Kallenberg |2002 , the diffusion process Z has a stationary distri- 
bution concentrated on M++ if and only if the process s(Z) has (— oo, +co) as its state 
space and the speed measure has finite total mass or s(Z) has a finite interval as its state 
space and the boundaries are reflecting. The introduction of an extra negative drift to 
geometric Brownian motion cannot make zero a reflecting boundary, so we are interested 
in conditions under which s(M++) = (—00,00) and the speed measure has finite total 
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mass. We see from (6.4) and (6.5) that this happens if and only if /j — a 2 /2 > 0, a 
condition we assume holds for the remainder of the proof. 

The diffusion s(Z) has a stationary distribution with density / := m (^^ + ^2 on 
s(M++) = (— oo,+oo), and so the stationary distribution of Z is the distribution on 
that has density 

g(x) = f(s(x))s'(x) 



1 



m(R ++ )a 2 (s(x)) 



s'(x) 



m(R ++ )a 2 (x)s'(x) 



m(R 



2 ^ 2 (f ) 



-2p/o* e |}(x-a) 



x £ 



This has the form of a Gamma(/c, 0) density with parameters 
Therefore, 



2k 



and k = ^ - 1. 



g{x) 



l 



r(fc)e j 



A-— 1 — £ 

-x e e 



1 



2/j. _ 

T 2 \ 



2k 



2/j 



e ^ 



Theorem 20.21 from Kallenberg 2002 implies that the shift-invariant cr-field is trivial 



for all starting points. The ergodic theorem for stationary stochastic processes then tells 
us that, if we start Z with its stationary distribution, 

2 nt roo 

lim - / h(Z s ) ds = h{x)g(x) dx 

t^co t Jq Jq 

for any Borel function h : R++ — > R with J* 0 °° \h(x)\g(x) dx < oo. Since Z has positive 
continuous transition densities we can conclude that 

ft roc 

h(Z s ) ds = h(x)g(x) dx 



lim - 

t— >oo t 



M -almost surely for any x £ 
In particular, 



xg(x) dx = k0 



ti- 



er 
~2 



Appendix B: Proof of Theorem 14.11 

To simplify our presentation, we re-write the joint dynamics of X and Y as 

(6.6) dX t = X t (fi-a-(aX t + cY t )) dt + a x X t dU t 

dY t = Y t (fi- (3 - (cX t + bY t )) dt + a Y Y t dV t , 

where a := (a, a) K , b := (/?, (3} K , c := (a, (3} K , ax '■= V& ■ Sa, and ay '■= VP ■ 



To prove Theorem 4.1 we need several preliminary results. First, we prove existence 
s of solutions to the sys 
Second, in Proposition 



6.1 



6.3 



Theorem 

(0, oo) 2 forall t > 0 whenever (Xq, Y 0 )e 



and uniq uene ss of solutions to the sy stem (6.6) as well as a useful comparison result in 

we establish that (Xt,Yt) remains in R\ + = 
Third, in Proposition 6.4, we show that 
weak limit points of the empirical measures | J Q * ¥^ x,y ^{(X s , Y s ) £ •} ds are stationary 
distributions for the process (X, Y) thought of as a process on R 2 ^ (rather than 
Finally, we show that lim^oo Yt = 0 with probability one in Proposition|6.5|and conclude 
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by showing that \ ^ X ^{{X S ,Y S ) G •} ds converges weakly to p x 
x {0}. 



5q concentrated 



on 



^4H 



The stochastic differential equation in (6.6) has a unique strong solution 



Theorem 6.1. 

and X t ,Y t e LP (P^)) for all t,p>0 for all (x, y) G M+ + . This solution satisfies X t > 0 
and Y t > 0 for all t > 0, pO"*) -almost surely for all (x,y) G R 2 ++ . Let ((X t ,Y t ))t>o be 
the stochastic process defined by the pair of stochastic differential equations 



(6.7) 



dXt = X t \ji ■ a — aXt) dt + a x X t dUt 
dY t = Y t (p-(3- bY t ) dt + oyYt dV t 



Lf(X 0 ,Y 0 ) = (X 0 ,Y 0 ), then 



and 



x t <x t 



Y t <% 



for all t > 0. 

Proof The uniqueness and existence of strong solutions is fairly standard, see, for ex- 
ample, Theorem 2.1 in Li and Mao |2009 . One notes that the drift coefficients are 
locally Lipschitz so strong solutions exist and are unique up to the explosion time. It 



is easy to show this explosion time is almost surely infinite (see Theorem 2.1 in Li and 
2009 ). Next, suppose that Xq = Xq. We adapt the comparison principle of Ikeda 



Mao 



and Watanabe (Chapter VI Theorem 1.1 from Ikeda and Watanabe 1989| ) proved by 



V.43.1 in |Rogers and Williams 2000 
Define p : R+ -)• R+ by p(x) 



the local time techniques of Le Gall ( see Theorem 1.4 from Le Gall [1983] and Theorem 



to show that Xt 
Note that 



X f > 0 for all t > 0. 



p{\x s - XsiyHiXs - x s > o} d[x - x] 



p(\X s - X 8 \)-\a x X B - a x X s ) 2 t{X s - X s > 0} ds 



<o 2 x t. 



Since L, p{u) 1 du = oo, by Proposition V.39.3 from Rogers and Williams 



2000 



the 

local time at 0 of X — X is zero for all t > 0. Put x + := x V 0. By Tanaka's formula 
(see equation IV.43.6 in Rogers and Williams 2000] ), 



(Xt-X t ) + = / t{X s - X s > 0}(a x X s - a x X s ) dU t 
Jo 



+ I t{X s - X s > 0} [(fi ■ a - (aX s + cY s ))X s - (p ■ a - aX s )X s ] ds. 
Jo 



For K > 0 define the stopping time 



T K := w£{t > 0 : X t > K or X t > K } 
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and the stopped processes X* = Xt k m and X L t 
at Tk and taking expectations yields 



K 



Xt k m- Then, stopping the processes 



0 < E(Xf - X 



E 



E 



rtAT K 

/ 1{X S - X s > 0} [(// • aX s - X s (aX s + cY s )) - (li ■ aX s - aX 2 s 
Jo 



ds 



tAT K 



< E 



t/\T K 



</x-aE 



1{X 
1{X 

t/\T K 



X s > 0} [li ■ a(X s - X s ) - a{X 2 s - X 2 S ) - cX s Y s ] ds 
X s > 0}li ■ a(X s - X s ) ds 
(X s -X s ) + ds 



</x-aE / (Xf -X 



ds. 



By Gronwall's Lemma (see, for example, Appendix 5 of Ethier and Kurtz [2005] ) 
E[(Xf - Xf )+]_= 0 for all t > 0, so Xf < Xf for all t > 0. Now let K -+ oo 
and recall that X does not explode to get that Xt < Xj for all i > 0. Since we have 
shown before that X is dominated by a geometric Brownian motion, a process that has 
finite moments of all orders, we get that Xt,Yt E L p (¥^ x ' y ' ) ) for all t,p > 0 and for all 
(x,y)£R 2 ++ . □ 

Remark 6.2. Note that the SDEs for all the processes considered here have unique strong 
solutions in LP for all t > 0, p > 0 and for all strictly positive starting points. This follows 
by arguments similar to those that are in Theorem 2.1 from Li and Mao |2009 and in 



Theorem 6.1 by noting that our SDEs for (X, Y), (X, Y) etc. are all of the form 



dX t 

dY t 

X 0 
% 



X t 
Yt 

x 

y 



Ai 
A 4 - 



- MY t 
^Xt 



A 3 Xt 
X 6 Y t 



dt + X t a x dU t 
dt + Y t a Y dV t 



for Ai, . . . , A6 € M + and x, y G 

The next proposition tells us that none of our processes hit zero in finite time. 



+> 



Proposition 6.3. Let (X,Y) be the process given by (6.6). Lf (Xq,Yq) G 
(Xt, Yt) G for all t > 0 almost surely. A similar conclusion holds for all of the 
other processes we work with. 

Proof. As an example of the method of proof, we look at the process (X, Y) given by 
(6.6). Taking logarithms and using Ito's lemma, 



Therefore, 



dlogX t 



logX t 



ix - a - (aX t + cY t ) - -o\ I dt + ax dU t . 



li - a - (aX s + cY s ) - -ax I ds + axU t - 



can't go to — oo in finite time because Xt and Yt do not blow up. 



□ 
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Proposition 6.4. Let (X, Y) be the process given by ( |6.6| ) and fix (x,y) G ^++- Any 
sequence {in}nsN such that t n — > oo has a subsequence {n n } ng pj such that the sequence 
of probability measures 



1 

u n JQ 



{(X S ,Y S ) G -}ds 



converges in the topology of weak convergence of probability measures on Any such 
limit is a stationary distribution for the process (X, Y) thought of as a process with state 
space 

Proof. Set ip(x,y) := x + y so that ip > 0 for x, y > 0. Put ip(x,y) = fj, • ax + fj, • 
(3y — x(ax + cy) — y(cx + by). Note that ip is bounded above on the quadrant x, y > 0 



and lim i 



(x,y)\\-too 



ip(x,y) = — oo where || • || is the Euclidean distance on M?. Using Ito's 



lemma we get 



ip{X s ,Y s )ds 



a Y Y s dV s + / a x X s dU s . 



Therefore, ip(Xt,Yt) — L ip(X s ,Y s ) ds is a martingale. Applying Theorem 9.9 of 



and Kurtz 1 2005 1 completes the proof. 



Ethier 



The following result is essentially Theorem 10 in Liu et al. [201 1 . We include the 
proof for completeness. 

Proposition 6.5. Suppose that a- /i— a-Sa/2 > 0, /3-fi— /3-S/3/2 > 0, and I (a, /3) < 0. If 

(X, Y) is the process given by (6.6), then lim^oo Yt = 0 ¥^ x ' y '-a.s. for all (x, y) G 

Proof. Using Ito's lemma and the definition of I(a,/3), 



log 



log 



a 



c fi • a 



Yt U t 
+ aa Y — - C(T x — 

al(a,f3)-(ab-c 2 )ti YsdS 



x 



+ aery 



V t 



u, 



ca x - 



By the Cauchy-Schwarz inequality, (ab — c 

+ l(a,f3) + a Y 



t t t 

2\ _ /„. „,\ ia a\ a„, a\ ^2 



t 



< 



{(3, /3) K — ((a, f5) K ) > 0, and so 
V t 



t 



t 



c U t 

-ox — - 
a t 



(6.8) 



Let X be the process defined by (6.7) with Xq = Xq. Proposition 

t 

2 



2.3 



implies 



lim - 

t— >co t 



X s ds = (fj, ■ a — o~ x /2) /a almost surely. 



o 



It follows from Theorem 6.1 that Xt < Xt for all t > 0. Thus, with probability one, 



logXj logA^ 
lim sup < lim sup 

£— >oo t ^on t 



t-, 


oo 




■ a 


-(> 


■ a 


= 0. 





a 



x 



1 f - Ut 
a lim - / X s ds + ax lim 



t— >oo t 



t— >oo t 



a 



x 



a(n ■ a — a x /2)/a 
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Since U and V are Brownian motions, lim^oo ^ = linn._s.oo ^ = 0, and 
lim sup^oo log t Xt < 0 almost surely, so 

log Yt 

lim sup < I(a, (3) < 0 almost surely. 

t— >oo * 

In particular, Hindoo Y t = 0 almost surely. □ 



We can now finish the proof of Theorem |4.1| Fix e > 0 and t] > 0 sufficiently small. 
Define the stopping time 

T e := inf{t > 0 : Y t > e}. 
and the stopped process X\ := X t /\T t - By Proposition 6.5, there exists T > 0 such that 

F (x,y){ Yt < e for all t > T} > 1 - rj 

Define the process X via 

dX t = Xt[(/i ■ a — ce) — aX t ] dt + oxX t dUt 

and the stopped process X\ := X t /\T e - Start the process X at time T with the condition 
Xt = Xt- We want to show that the process X e is dominated by the process X e , that 
is X\ > XI for all t >T. By the strong Markov property, we can assume T = 0. 

The proof is very similar to the one from Theorem 6.1 With the notation from the 
proof of Theorem 6.1 we have 

/ P (\xi - xi\)-H{xi - xt > o} d[x* - x% 

Jo 

= [ p(\Xl - Xl\)~\a x Xt - a x Xl)H{Xl - XI > 0}] ds 
Jo 

<a\t 

so the local time of the process X e — X e at zero is identically zero. Then, using Tanaka's 
formula, 

rtAT e 

(X e t -X f t )+ = / t{X s -X s >0}{a x X s -a x X s )dU t 
Jo 

rtAT e 

+ / t{X„ - X s > 0} [((// • a - ce)X s - aX 2 s ) - (fi ■ aX s - X s (cY s + aX s ))] ds. 
Jo 

Taking expectations, 

rtAT e 

E[(X e t - Xl)+] = E t{X s -X s >0}[(p-a(X 3 -X s )-{ceX s -cX s Y s ) 
Jo 

-a(X 2 s -X 2 s ))ds] 
rtATe 

< p-aE (X s -X s ) + ds 
Jo 



< fi-aE {XI -X e s ) + ds. 



o 



By Gronwall's Lemma, IE[(-_| — X£y] =0. As a result, remembering we assumed 
T = 0, we have X\ < Xf for all t > T. For e small enough we know that X has a 
stationary distribution concentrated on For any sequence a n — > oo, if the Cesaro 

averages ^- J^ n F^ x ' y ^{(X s , Y s ) _ •} ds converge weakly, then the limit is a distribu- 
tion of the form ip (8) 5q, where p> is a mixture of the unique stationary distribution 
p x described in Proposition 2.3 and the point mass at 0. By the above, the limit of 
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T~ Jo™ ^ {(X s ,Y s ) G -}ds cannot have any mass at (0,0) because X t < X t on the 
event {Y t < e for all t > T} that has probability F^{Y t < e for all t > T} > 1 - rj. 
Since rj > 0 was arbitrary, we conclude that <p = p x <8> So, as required. 



Appendix C: Proof of Theorem 14.21 



Our proof is along the same lines as the proofs of Theorems 4 and 5 in Schreiber et al. 



2011 



We will once again simplify our notation by re-writing the SDE for the pair (X, Y) 
as in (6.6). We assume throughout this appendix that the hypotheses of Theorem |4.2| 
hold; that_is,_X(a, (3) > 0 and Z(/3, a) > 0. 

Let ((Xt,Yt))t>o be the stochastic process defined by the pair of stochastic differ- 
ential equations in (6.7) with initial conditions (Xq,Yq) = (Xq,Yq). We know from 
Theorem [6t1 that X t < X t and Y t < Y t for all t > 0. 

Note from Corollary 3.3 that a ■ (fx — Sa/2) > 0 and /3 ■ (fi — S/3/2) > 0 and hence, 
by Proposition 2.3 the process (X,Y) has a unique stationary distribution on and 
is strongly ergodic. 
Let 



t{(X s ,Y s ) e -}ds 



be the normalized occupation measures of (X, Y). We know that the random probability 
measures 

1 



n*(0 := 



t 



i{(x s ,y s ) e -}ds 



o 



converge almost surely and so, in particular, they are tight on 
any e > 0 we can find a box [Q,K] x [0, K] such that 

ft 

±{(X S , Y s ) £ [0, K] x [0, K]} ds > 1 - e for all t > 0. 



0, oo) 2 ; that is, for 



Therefore, 



1{(X S ,Y S ) e [0,K] x [0,K]}ds > 1{(X S ,Y S ) G [0,K] x [0,K]}ds 

t Jo 

> 1 - e for all t > 0, 



and hence the normalized occupation measures of (X, Y) are also tight on 



By 



Prohorov's theorem Kallenberg, 2002, Theorem 16.3], there exists a random probability 
measure v on and a (possibly random) sequence (t n ) C M ++ such that t n — > oo for 
which 



(6.9) 



n 



V 



as n — > oo almost surely, where ==> denotes weak convergence of probability measures on 
M 2 ^. That is, with probability one for all bounded and continuous function u : — >■ M 
we have 



u(x,y) U tn (dx,dy) 



u(x,y)v(dx,dy) 



as n 



oo. 



Proposition 6.6. The probability measure v is almost surely a stationary distribution 
for (X,Y) thought of as a process with state space M 2 ^. 
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Proof. Let (Pt)t>o be the semigroup of the process (X, Y) thought of as a process on 
M+. For simplicity let us write Z t := (X t , Y t ) for all t > 0 and v n := Ut n . 

By the Strong Law of Large Numbers for martingales, we have that for all r G R + 
and all bounded measurable functions / 

j fe-i 

} im t. ^2if( Z r+{i+i)t) ~ Ptf{Z r +u)] = 0 almost surely. 



As a result, 



1 

k 



kt 



f(Z s ) ds 



(k-i)t 



k-l 



Ptf(Z s )ds 



i=0 



/ [f(Z r+ (i + i)t) - Ptf{Z r+it )\ dr 



0 as k — > oo almost surely. 



This implies that 



lim — 

u— >oo u 



pu 

/ [f(Z s+t ) - Pt{Z s )} ds = 0 almost surely. 

JO 



Thus, 



Ptfdv 



lim 

n— >oo 



ffdu n -f 



Ptfdv, 



lim — 
lim — 



(f(Z s )-P t f(Z s ))ds 



tn t 



(f(Z s+t ) - P t f{Z s )) ds + / f{Z s ) ds 



Ptf(Z s ) ds 



lim 

n— >co t 



1 



t n -t 



(6.10) 



(f(Z s+t )-P t f(Z s ))ds 
0 almost surely. 



The last result is equivalent to saying that v is almost surely a stationary distribution 
for (X,Y). 

□ 

Proposition 6.7. There exists a stationary distribution ir of (X,Y) that assigns all of 



its mass to 



Proof. We argue by contradiction. Because the process stays in one of the four sets 

x {0}, {0} x {(0,0)} when it is started in the set, any stationary distribution 

for (X, Y) thought of as a process on can be written as a convex combination of 
stationary distributions that respectively assign all of their masses to one of the four 
sets, should such a stationary distribution exist for the given set. Suppose there is no 
stationary distribution that is concentrated on Then, any stationary distribution 

is the convex combination of stationary distributions that respectively assign all of their 
mass to the three sets M++ x {0}, {0} x and {(0,0)}, and hence any stationary 

distribution is of the form 

PXVX+PYVY +R>Vo)> 
where the random variables px ,Py,Po are nonnegative and px + Py + Po = 1 almost 
surely, and fix = Px ® and Hy = ® Py for p x and p Y the unique stationary 
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distributions of X and Y. Next, we proceed as in Proposition 6.5 to find the limit of 
lQ g^t n ^ ug £ rg ^. ar g ue that 



(6.11) 



lim 

n— >oo t 



— f " X s ds = [ 

tn Jo Jm. 2 , 



lim 

n— >oo t 



1 f tn f 

— / Y s ds = I yv(dx,dy) almost surely. 
tn Jo Jr2 



Note that the infinitesimal generator of (log X, log Y) thought of as a process on M? 
is uniformly elliptic with smooth coefficients and so it has smooth transition densities 
(see, for example, Section 3.3.4 of Stroock |2008|). Moreover, an application of a suitable 
minimum principle for the Kolmogorov forward equation (see, for example, Theorem 5 
in Section 2 of Chapter 2 of Friedman [1964] ) shows that the transition densities are 
everywhere strictly positive. It follows that (X, Y) thought of as a process on has 
smooth transition densities that are everywhere positive. 

Because the process X also has smooth, every positive transition densities for similar 
reasons, the almost sure behavior of the X started from a fixed point is the same as it is 
starting from its stationary distribution p%- . As a result, we get by Birkhoff's pointwise 
ergodic theorem Kallenberg, 2002, Theorem 10.6] that, for all K > 0, 



lim 



1 



n— >oo t 



X s t{X s > K}ds = E p * [X s l{X s > K}] 



n JO 



almost surely for any x G M+. Therefore, by dominated convergence 



lim lim — 

K^roo n— too t n Jq 



X s t{X s > K}ds= lim E p JX s t{X s > K}] = 0. 



The following inequalities are immediate due to the positivity of the terms 



1 

tn 



X s t{X s <K}ds< 



3.12) 



1 

tn Jo 
1 

tn Jo 
1 



+ 



t 



n JO 



Recall that Xt < Xt for all t > 0 and hence 
1 

tn Jo 



tn 1 

X S 1{X S >K}ds<- 



Xe ds 



X S 1{X S < K} ds 



X S 1{X S > K}ds. 



X s t{X s > K}ds. 



n JO 



This implies 



lim sup — 

n— >oo t n 



1 



X S 1{X S > K}ds < lim sup 

n— >oo t n Jo 



X S 1{X S > K}ds, 



and therefore 



3.13) 



0 < lim lim sup — 

K— >co n — j^oo t n Jq 

1 

< lim lim sup — 

K— >co n — j^oo t n Jq 



X s l{X s > K}ds 
X s t{X s > K}ds = 0. 
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By ([679]) and Theorem 4.27 of |Kallenberg| (2002 



lim — f " X s t{X s < K} ds 

n->oo t n Jq 



xt{x < K} u{dx, dy). 



for any K such that 

v({K} x R + ) = 0. 

While this last condition need not hold a priori for all K, we can only have 

v{{K} x M + ) > 0 

for countably many K, so there exists a sequence (K m ) C M+ such that K n 
m — > oo with 

u{{K m } x M+) = 0. 

By dominated convergence, 



oo as 



(6.14) 



lim lim — 



X S 1{X S < K m } ds = lim / xl{x < K} v{dx, dy) 

K^oo Its/2 



x u(dx, dy). 



Combining (6.12), (6.13) and (6.14) gives (6.11). 



It follows from Ito's formula, the observation I(a,a) = 0, (6.11), and the fact that 



lim 



u t , 

n— >oo + 



-- 0 that 
loeXi 



lim 

n— »oo t r 



fJi . a -^-E v [aX t + bY t ] 

2 

Px { u-a-aW^iXt)-^- 
+ py {^■a-b-E^[Y t ]- a f 



+ Po M " a 



a 



x 



PxZ(a, a) + pyZ(a, (3) + p 0 [ p ■ a - 



pyT{a, (3) + po [ n ■ a - 



almost surely. 



By assumption, X(a, 0) > 0 and we have already observed that [i-a — ^if- > 0. Because 
Xt converges in distribution as t — > oo to a distribution that assigns all of its mass to 
it follows that log ^ tn converges in probability to 0. However, since Xt < Xt for 

all t > 0 it follows that pyX(a, f3) + po ^/i • a — < 0 and hence 

(6.15) py = Po = 0 almost surely. 
The same argument applied to (Y t )t>o establishes 

(6.16) px = Po = 0 almost surely. 

Therefore, px = Py = Po = 0, and this contradicts the assumption that px +Py +Po = 
1. □ 



We can now finish the proof of Theorem 4.2 
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Proof. Proposition 6.7 implies that (X,Y) has a stationary distribution tt on By 
Theorem 20.17 from Kallenberg |2002] , our process (X,Y) is either Harris recurrent or 
uniformly transient. We say that (X t , Yt) — ¥ oo almost surely as t — > oo if tx j Xt, Y t ) — > 0 
oo for any compact set K C M.2 i_ ■ Theorem 20.21 from Kallenberg 



2002 



gives 



oo and so (X, Y) cannot have a stationary 



as t 

that if (X,Y) is transient, then (Xt,Yt 
distribution. Hence, since we know our process has a stationary distribution tt, it must 
be Harris recurrent. Theorem 20.21 from Kallenberg 12002] then gives us equation (4.1 ). 



Theorem 20.18 from Kallenberg 2002 , 20.18 gives that any Harris recurrent Feller 
process on with strictly positive transition densities has a locally finite invariant 
measure that is equivalent to Lebesgue measure and is unique up to a normalization. We 
already know that we have a stationary distribution, so this distribution is unique and 
has an almost everywhere strictly positive density with respect to Lebesgue measure. 
Theorem 20.12 from |Kallenberg| [2002| says that any Harris recurrent Feller process is 
strongly ergodic, and so equation (4.2) holds. 

□ 



Remark 6.8. In Theorem 3.1 of Zhang and Chen 20131, the authors claim to show that 



the system of SDE describing (X, Y) always has a unique stationary distribution. We 
note that their use of moments just checks tightness in := [0, oo) 2 and not in R\+ = 
(0, oo) 2 . It does not stop mass going off to R 2 + \ R% + = (R + x {0}) U ({0} x R + ), which 
is exactly what can happen in our case. Thus, their proof only shows the existence of a 
stationary distribution on - it does not show the existence of a stationary distribution 
on Furthermore, their proof for the uniqueness of a stationary distribution on 

breaks down because their assumption of irreducibility is false. The process (X, Y) is 
irreducible on but it is not irreducible on M 2 . since P t ((0, 0), U) := P(°'°){(X t , Y t ) £ 
U} = 0 for any open subset U that lies in the interior of If we work on it is 
not true that the diffusion (X, Y) has a unique stationary distribution. We can obtain 
infinitely many stationary distributions on R 2 ^ of the form (up x + v5q) £g> 5q where p x 
is the unique stationary distribution of X on M++ and u,v £ M + satisfy u + v = 1. 

Appendix D: Proof of Theorem 15.11 

Assume that the matrix S is positive definite and that the dispersion proportion 
vector a is such that p, ■ a — a ■ T,a/2 > 0 so that a population playing the strategy a 
persists. Under these assumptions the function f3 h-> I(a, f5) is strictly concave. Hence, 
by the method of Lagrange multipliers, I(a, j3) < 0 for all f3 ^ a and an > 0 for all i if 
and only if there exists a constant, which we denote by A, such that 

dl 



3.17) 



A 



0A 



/3=a 



Hi - KiCXi(n ■ a- a - Sa/2)/ (a, a) K - 2^ ajOij 



for all i. Multiplying (6.17) by ctj and summing with respect to i, we get 

Ea/2)/(a,a) K 



A 



fx - a — (a, a) K (ji ■ a 
—a ■ Sq/2 



a 



a ■ Sa 



This expression for the Lagrange multiplier and (6.17) provide the characterization of 



a mixed ESS in equation (5.1) when a« > 0 for all i. The characterization of the more 



general case of ai > 0 for at least two patches follows similarly by restricting the method 
of Lagrange multiples to the appropriate face of the probability simplex. 

Suppose that pi — au/2 > 0 so that a population remaining in patch i and not 
dispersing to other patches persists. The strategy aj = 1 and a,- = 0 for all j ^ i is an 
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ESS only if 

< 0 



for all j 7^ i. Evaluating these partial derivatives gives the criterion (5.2) for the pure 
ESS. 

We conclude by considering the case n = 2. Define the function g : [0, 1] — )■ M by 



dl 

= ^-(( a l) a 2) 5 (&1,&2)) 

opi 
dl 

~ 7M"(( a l' a 2)' ( & 1' & 2)) 
CP2 



(ai,a2)=(a,l— a),(bi,6 2 )=(a,l-a) 



(ai,a 2 )=(a,l— a),(6i,6 2 )=(a,l-a) 



The inequalities (5.2) for the pure strategies (1,0) and (0,1), respectively, correspond 
to g(0) < 0 and g(l) > 0, respectively. Hence, when these inequalities are reversed, the 
intermediate value theorem implies there exists a E (0, 1) such that g(a) = 0. Such an 



a satisfies the mixed ESS criterion (5.1) and, therefore, is an ESS. 
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